Universidad de Las Palmas de Gran Canaria
Zaragoza 3 de diciembre de 2019
Grandes retos de la utilización de datos masivos
1. Volumen de información
2. Velocidad de generación de datos
3. Variedad de datos
4. Veracidad de los datos
5. Viabilidad del uso de los datos
6. Visualización de los datos
7. Valor de los datos: datos que se transforman en información, conocimiento, acción y decisión
Grandes retos de la utilización de datos masivos
1. Volumen de información4. Veracidad de los datos
5. Visualización de los datos
6. Viabilidad del uso de los datos
7. Valor de los datos: datos que se transforman en información, conocimiento, acción y decisión
Presentar las oportunidades que ofrece el CMBD en la investigación en servicios sanitarios.
Efecto de la adopción de trombectomía primaria en la tasa de mortalidad por accidente cerebrovascular agudo (Indicador de calidad IQI 17).
Datos individuales de ingresos hospitalarios ente Enero de 2003 y diciembre de 2015.
Pacientes que ingresan con diagnóstico principal ictus isquémico
“Un término ambiguo y potencialmente engañoso que se refiere a los datos recopilados regularmente (a veces de manera rutinaria, pero a menudo con un esfuerzo y gasto considerable para maximizar la VALIDEZ), principalmente por organizaciones gubernamentales, este es el caso de datos sobre nacimientos, defunciones, registros de enfermedades o estadísticas de servicios sanitarios.”
- Hospitalización (CMBD-H)
Información administrativa y clínica de los pacientes ingresados en unidades de hospitalización, la cual es recogida y codificada sistemáticamente en los hospitales, para cada episodio de hospitalización.
Los registros de alta son enviados al Ministerio por las comunidades autónomas, de acuerdo con el modelo de datos establecido al efecto por el Consejo Interterritorial.
Datos normalizados desde 1997 hasta 2015 inclusive de todos los hospitales de agudos del Sistema Nacional de Salud.
Las variables de consulta incluyen las básicas sobre edad y sexo del paciente, del episodio (tipo de ingreso y alta) y variables clínicas (sobre diagnósticos y procedimientos) a las que se añaden otras derivadas del uso de sistemas de clasificación de pacientes.
Flowchart of the patient population meeting inclusion criteria for the study
| Variable | Year.2003 | Year.2006 | Year.2009 | Year.2012 | Year.2015 | p_value |
|---|---|---|---|---|---|---|
| Number of admissions | 12,976 | 13,505 | 15,227 | 16,829 | 17,658 | <0.001 |
| In-hospital mortality | 1,686 | 1,647 | 1,735 | 1,809 | 1,848 | <0.001 |
| Age, year (SD) | 73.59 (12.17) | 73.46 (12.75) | 73.68 (13.04) | 74.27 (13.11) | 74.04 (13.44) | <0.001 |
| Gender, female (%) | 46.59 | 46.60 | 46.95 | 47.59 | 46.15 | 0.252 |
| Public holiday admission (%) | 26.88 | 27.09 | 26.98 | 27.71 | 27.44 | 0.001 |
| Congestive heart failure (%) | 4.09 | 4.15 | 4.55 | 5.41 | 5.41 | <0.001 |
| Pulmonary circulation disease (%) | 1.03 | 1.21 | 1.60 | 1.77 | 2.07 | <0.001 |
| Chronic pulmonary disease (%) | 9.65 | 8.34 | 8.37 | 9.30 | 9.44 | 0.741 |
| Renal failure (%) | 2.89 | 3.47 | 5.62 | 6.83 | 7.89 | <0.001 |
| Lymphoma (%) | 0.22 | 0.33 | 0.35 | 0.30 | 0.42 | 0.015 |
| Metastatic cancer (%) | 0.94 | 1.20 | 1.18 | 1.36 | 1.44 | <0.001 |
| Coagulopathy (%) | 0.62 | 0.70 | 0.72 | 1.31 | 1.26 | <0.001 |
| Fluid and electrolyte disorders (%) | 1.52 | 1.47 | 1.90 | 2.35 | 2.83 | <0.001 |
| Drug abuse (%) | 0.28 | 0.59 | 0.72 | 0.66 | 0.83 | <0.001 |
| Endovascular treatment (%) | 0.52 | 2.49 | 7.00 | 12.10 | 15.45 | <0.001 |
Tenemos una variable respuesta y (Mortalidad intrahospitalaria por ictus isquémico) la cual queremos estudiar. Para ello seleccionamos una muestra de tamaño n de nuestra población de interés. Deseamos inferir propiedades de la variable y en términos de k predictores observados \[X_i (x_{1i},.....x_{ki})\] \[i=1,....,n\] Usamos estas k variables para explicar la variación en y
Nuestra variable endógena Y representa la Mortalidad en el hospital y sólo puede tomar dos valores, que codificamos como \(1\)= alta exitus y 0 para el resto. La distribución de Y es una Bernoulli cuya esperanza es:
\[E[Y]=P[Y=1]=p (0<p<1)\]
Sea X una posible variable predictora de Y, entonces la distribución condicional de Y sobre un valor de X = x, también sigue una distribución de Bernoulli de forma que la esperanza condicionada de Y sobre X = x es:
\[ E[Y |X = x] = P [Y = 1 |X = x] = p(x) \] y la varianza condicionada
\[ Var[Y |X = x] = p(x)·(1-p(x)) \]
Podemos formular el modelo para Y en función de X si para su estimación utilizamos una regresión lineal clásica. El modelo quedaría como: \[ Y = \alpha + \beta x + \epsilon(x) \] donde los errores son variables aleatorias independientes con esperanza 0 que distribuyen con una Bernoulli.
Esta estimación se corresponde con el denominado modelo de probabilidad lineal.
##
## Call:
## lm(formula = Yvar ~ Age + Pub.Hol + chf + pulmcirc + chrnlung +
## renlfail + lymph + mets + coag + lytes + procedure + Trend +
## as.factor(Season), data = df.bis)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.64656 -0.13537 -0.10047 -0.04178 1.12333
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.589e-01 4.804e-03 -33.080 < 2e-16 ***
## Age 3.900e-03 5.441e-05 71.674 < 2e-16 ***
## Pub.Hol1 4.326e-03 1.562e-03 2.769 0.005616 **
## chf1 1.363e-01 3.326e-03 40.970 < 2e-16 ***
## pulmcirc1 6.846e-02 5.621e-03 12.180 < 2e-16 ***
## chrnlung1 1.060e-02 2.427e-03 4.368 1.26e-05 ***
## renlfail1 3.743e-02 3.085e-03 12.133 < 2e-16 ***
## lymph1 2.741e-02 1.220e-02 2.247 0.024663 *
## mets1 2.148e-01 6.240e-03 34.417 < 2e-16 ***
## coag1 4.256e-02 7.392e-03 5.757 8.56e-09 ***
## lytes1 1.502e-01 5.071e-03 29.620 < 2e-16 ***
## procedure 2.317e-02 2.662e-03 8.706 < 2e-16 ***
## Trend -2.818e-04 1.594e-05 -17.685 < 2e-16 ***
## as.factor(Season)2 -8.379e-03 3.403e-03 -2.462 0.013820 *
## as.factor(Season)3 -6.544e-03 3.326e-03 -1.968 0.049116 *
## as.factor(Season)4 -1.152e-02 3.373e-03 -3.415 0.000639 ***
## as.factor(Season)5 -1.814e-02 3.344e-03 -5.426 5.78e-08 ***
## as.factor(Season)6 -1.614e-02 3.365e-03 -4.795 1.63e-06 ***
## as.factor(Season)7 -2.141e-02 3.367e-03 -6.360 2.02e-10 ***
## as.factor(Season)8 -1.456e-02 3.405e-03 -4.276 1.90e-05 ***
## as.factor(Season)9 -1.899e-02 3.415e-03 -5.562 2.67e-08 ***
## as.factor(Season)10 -1.226e-02 3.353e-03 -3.657 0.000256 ***
## as.factor(Season)11 -1.288e-02 3.346e-03 -3.850 0.000118 ***
## as.factor(Season)12 -5.443e-03 3.315e-03 -1.642 0.100567
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3078 on 196075 degrees of freedom
## Multiple R-squared: 0.0542, Adjusted R-squared: 0.05409
## F-statistic: 488.5 on 23 and 196075 DF, p-value: < 2.2e-16
| Variables | Coeff.est | Variables.continuacion | Coeff.esti |
|---|---|---|---|
| Age | 0.15982859 | Coagulopathy | 0.01266760 |
| Congestive heart failure | 0.09132888 | Aug | -0.01241610 |
| Metastatic cancer | 0.07571767 | Nov | -0.01133821 |
| Fluid and electrolyte disorders | 0.06540839 | Oct | -0.01074906 |
| Trend | -0.03985476 | Apr | -0.00997658 |
| Renal failure | 0.02692666 | Chronic pulmonary disease | 0.00961470 |
| Pulmonary circulation disease | 0.02691802 | Feb | -0.00714759 |
| Trombechtomic procedure | 0.01952882 | Public Holiday | 0.00608986 |
| Jul | -0.01862153 | Mar | -0.00581105 |
| Sep | -0.01612665 | Lymphoma | 0.00493556 |
| May | -0.01596039 | Dec | -0.00486772 |
| Jun | -0.01404541 |
Los GLM’s tienen tres componentes.
Los modelos logit son modelos GLM dónde la distribución de la variable \(Y\) es binomial y la función de enlace es la función logística. Si llamamos \(\mu=P(Y=1)\) el modelo queda especificado como:
\[logit[\mu]=log\left(\frac{\mu}{1-\mu}\right)=\alpha+\beta_1x_1+...+\beta_px_p\] Estimemos un modelo logit sobre nuestros datos…
##
## Call:
## glm(formula = Yvar ~ Age + Pub.Hol + chf + pulmcirc + chrnlung +
## renlfail + lymph + mets + coag + lytes + procedure + Trend +
## as.factor(Season), family = binomial, data = df.bis)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.9269 -0.5274 -0.4246 -0.2882 3.2917
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -6.0112814 0.0648253 -92.730 < 2e-16 ***
## Age 0.0537276 0.0007529 71.363 < 2e-16 ***
## Pub.Hol1 0.0436500 0.0163794 2.665 0.007700 **
## chf1 0.8596751 0.0256122 33.565 < 2e-16 ***
## pulmcirc1 0.5204587 0.0461909 11.268 < 2e-16 ***
## chrnlung1 0.1317634 0.0239292 5.506 3.66e-08 ***
## renlfail1 0.3143656 0.0277984 11.309 < 2e-16 ***
## lymph1 0.2874196 0.1143866 2.513 0.011981 *
## mets1 1.5852291 0.0462862 34.248 < 2e-16 ***
## coag1 0.4016866 0.0674724 5.953 2.63e-09 ***
## lytes1 0.9292229 0.0380178 24.442 < 2e-16 ***
## procedure 0.3033109 0.0281435 10.777 < 2e-16 ***
## Trend -0.0033268 0.0001698 -19.591 < 2e-16 ***
## as.factor(Season)2 -0.0795154 0.0347926 -2.285 0.022289 *
## as.factor(Season)3 -0.0587290 0.0339390 -1.730 0.083554 .
## as.factor(Season)4 -0.1089387 0.0347928 -3.131 0.001742 **
## as.factor(Season)5 -0.1804679 0.0349378 -5.165 2.40e-07 ***
## as.factor(Season)6 -0.1645463 0.0351628 -4.680 2.88e-06 ***
## as.factor(Season)7 -0.2191299 0.0356522 -6.146 7.93e-10 ***
## as.factor(Season)8 -0.1479730 0.0354570 -4.173 3.00e-05 ***
## as.factor(Season)9 -0.1986191 0.0360025 -5.517 3.45e-08 ***
## as.factor(Season)10 -0.1217769 0.0346372 -3.516 0.000438 ***
## as.factor(Season)11 -0.1269012 0.0346102 -3.667 0.000246 ***
## as.factor(Season)12 -0.0530982 0.0337046 -1.575 0.115164
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 138278 on 196098 degrees of freedom
## Residual deviance: 127665 on 196075 degrees of freedom
## AIC: 127713
##
## Number of Fisher Scoring iterations: 5
##
## Attaching package: 'MASS'
## The following object is masked from 'package:plotly':
##
## select
## Waiting for profiling to be done...
## 2.5 % 97.5 %
## (Intercept) 0.0025 0.0022 0.0028
## Age 1.0552 1.0536 1.0568
## Pub.Hol1 1.0446 1.0116 1.0786
## chf1 2.3624 2.2464 2.4837
## pulmcirc1 1.6828 1.5362 1.8412
## chrnlung1 1.1408 1.0884 1.1954
## renlfail1 1.3694 1.2965 1.4457
## lymph1 1.3330 1.0595 1.6598
## mets1 4.8804 4.4554 5.3419
## coag1 1.4943 1.3071 1.7030
## lytes1 2.5325 2.3500 2.7277
## procedure 1.3543 1.2813 1.4308
## Trend 0.9967 0.9963 0.9970
## as.factor(Season)2 0.9236 0.8626 0.9887
## as.factor(Season)3 0.9430 0.8823 1.0078
## as.factor(Season)4 0.8968 0.8376 0.9600
## as.factor(Season)5 0.8349 0.7796 0.8940
## as.factor(Season)6 0.8483 0.7917 0.9087
## as.factor(Season)7 0.8032 0.7489 0.8613
## as.factor(Season)8 0.8625 0.8045 0.9245
## as.factor(Season)9 0.8199 0.7639 0.8797
## as.factor(Season)10 0.8853 0.8272 0.9475
## as.factor(Season)11 0.8808 0.8230 0.9426
## as.factor(Season)12 0.9483 0.8876 1.0130
- AIC (Criterio de información de Akaike):
métrica análoga al \(R^2\) ajustado. Siempre preferimos el modelo con un valor mínimo de AIC.
- Null Deviance y Residual Deviance:
suma de los cuadrados de los residuos de la devianza en el modelo con solo constante y en el modelo con covariables respectivamente. Cuando menor valor mejor es el modelo.
- Valores observados vs. predichos:
porcentaje de observaciones correctamente predicas;
sensibilidad;
especificidad
área bajo la curva ROC.
## [1] "AIC 127713.4842"
## [1] "RsqrMcFadden 0.0767"
## [1] "Residual deviance 127665.4842"
##
## FALSE TRUE
## 0 0.549824 0.337258
## 1 0.035798 0.077119
## [1] "Accuracy 0.6269"
## [1] "Sensibilidad 0.683"
## [1] "Especificidad 0.6198"
## [1] "Area under the ROC curve"
## [1] 0.70728
## Analysis of Deviance Table
##
## Model: binomial, link: logit
##
## Response: Yvar
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 196098 138278
## Age 1 7123 196097 131155 < 2e-16 ***
## Pub.Hol 1 6 196096 131149 0.014 *
## chf 1 1214 196095 129935 < 2e-16 ***
## pulmcirc 1 129 196094 129806 < 2e-16 ***
## chrnlung 1 35 196093 129771 3.4e-09 ***
## renlfail 1 105 196092 129666 < 2e-16 ***
## lymph 1 6 196091 129660 0.012 *
## mets 1 948 196090 128712 < 2e-16 ***
## coag 1 32 196089 128680 1.5e-08 ***
## lytes 1 496 196088 128184 < 2e-16 ***
## procedure 1 46 196087 128138 1.3e-11 ***
## Trend 1 399 196086 127739 < 2e-16 ***
## as.factor(Season) 11 74 196075 127665 2.3e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
En los modelos de regresión es habitual asumir una relación lineal entre la variable dependiente y las variables independientes.
Como primera mejora se propone una relación polinómica, la cual obtiene mejores resultados la mayoría de las veces, pero que en situaciones de elevada variabilidad lleva a un ajuste excesivo que deriva en importantes errores de predicción.
Una alternativa a la relación polinómica es, en lugar de estimar un polinomio sobre todo el conjunto de valores de la variable (por ejemplo, la edad), dividir la misma en tramos ajustando relaciones diferentes lineales/polinómicas etc. en cada uno de ellos.
Charles J. Stone (1985). Additive Regression and Other Nonparametric Models.
The Annals of Statistics Vol. 13, No. 2, pp. 689-705.
Para predecir a partir de complejas relaciones (no lineales) que podrían incluso interaccionar;
Para comprender, visualizar esas relaciones;
Para controlar por dichas relaciones.
-Generalizado: pueden manejar muchas distribuciones de datos (normales, binomiales, de recuento, etc.);
-Aditivo: los términos simplemente se suman, si bien los términos en sí no son lineales;
-Modelo: Modelo
Un modelo de regresión no-lineal general tendría la siguiente forma:
\[ y=f(\mathbf{X}\beta) + \epsilon\] La función \(f(.)\) relaciona el valor medio de la variable respuesta \(y\) con las variables predictoras, y su forma está especificada a priori.
Los modelos GAM permiten otra generalización al incorporar formas no lineales diferentes para cada uno de los predictores.
\[ y=\beta_0+f_1(x_{i1})+f_2(x_{i2})+...+f_p(x_{ip})+\epsilon_i\] ajusta una función \(f_j\) por separado para cada \(X_j\). En otras palabras \(\beta x\) se reemplaza por la función de suavizado \(f(x)\).
Regla flexible de dibujo técnico
En los modelos GAM pasamos de los splines de dibujo técnico a los splines de regresión. Las regresiones splines (loess, etc.) no conectan todos los puntos sino que los ajustan a una función de suavizado simple en términos de grados de libertad.
El ajuste de los datos en cada sector se realiza de forma independiente.
La función de suavizado definitiva se elige teniendo en cuenta el “equilibrio” entre la mejoría en la verosimilitud (por ejemplo, la reducción en el AIC) y los grados de libertad.
###Modelo de regresión aditivo Si las \(X'_s\) fuesen completamente independientes podríamos estimar una regresión no-paramétrica para cada una de ellas por separado.
Pero como cabe esperar que las \(X'_s\) estén relacionadas necesitamos ir “eliminando” el efecto de cada predictor.
Se utiliza el algoritmo denominado backfitting.
En un modelo con dos variables: \[ y=\beta_0+f_1(x_{i1})+f_2(x_{i2})+\epsilon_i\] reagrupamos la ecuación como \[ y-f_2(x_{i2})=\beta_0+f_1(x_{i1})+\epsilon_i\] en otras palabras “conocer una función nos permite encontrar la siguiente”
Se recalculan las variables como desviaciones respecto de su media, eliminando de esta manera los interceptos \(\beta_0\).
Como estimaciones preliminares para cada función se utiliza la estimación mínimos cuadrados de \[y^*_i=b_1x^*_{i1}+b_2x^*_{i2}+\epsilon_i\]
A partir de las regresiones parciales obtengo los residuos parciales de cada \(X\).
En el ejemplo de dos covariables, los residuos parciales de \(X_1\) son \[e^{(1)}_{i[1]}=y^*_i-b_2(x^*_{i2})\] y lo mismo con los de \(X_2\)
A continuación, suavizamos estos residuos respecto de sus respectivas \(X\), obteniéndose nuevos suavizados de los residuos parciales corregido el efecto del resto de \(X\).
\[RSS(f,\lambda)=\sum_{i=1}^N (y_i-f(x_i))^2+\lambda\int_{x_1}^{x_n}f''(x)^2~dx\] el primer término mide la cercanía a los datos, y el segundo penaliza la curvatura en la función.
El proceso de búsqueda de nuevas estimaciones de las funciones de suavizado de los residuos parciales se reitera hasta la convergencia
Cuando termina el proceso anterior obtenemos estimaciones de \(s_j(X_{ij})\) para cada valor de \(X_j\)
Lo "más interesante, hemos reducido una regresión múltiple a una serie de regresiones parciales bidimensionales, facilitando la interpretación (5. Visualización). Dado que cada regresión parcial es solo bidimensional, las formas funcionales se pueden trazar en gráficos bidimensionales que muestran los efectos parciales de cada \(X_j\) en \(Y\).
Incorporamos las funciones de suavizados sobre nuestras variables continuas: Edad, tendencia y estacionalidad.
\[logit[\mu]=log\left(\frac{\mu}{1-\mu}\right)=\alpha+f_1(x_1)+\beta_2x_2+...+f_p(x_p)\]lm(y ~ x1 + x2, data=data)
glm(y ~ x1 + x2, data=data, family=binomial)
library(mgcv)
gam(y ~ x1 + s(x2), data=data, family=binomial , method = “REML”)
##
## Family: binomial
## Link function: logit
##
## Formula:
## Yvar ~ s(Age) + Pub.Hol + chf + pulmcirc + chrnlung + renlfail +
## lymph + mets + coag + lytes + procedure + s(T) + s(Nmonth,
## k = 12, bs = "cc")
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.4254 0.0108 -225.20 < 2e-16 ***
## Pub.Hol1 0.0450 0.0164 2.73 0.0062 **
## chf1 0.8407 0.0259 32.50 < 2e-16 ***
## pulmcirc1 0.5324 0.0465 11.46 < 2e-16 ***
## chrnlung1 0.1648 0.0241 6.84 8.0e-12 ***
## renlfail1 0.3207 0.0280 11.45 < 2e-16 ***
## lymph1 0.3313 0.1148 2.88 0.0039 **
## mets1 1.6142 0.0461 35.00 < 2e-16 ***
## coag1 0.4153 0.0676 6.14 8.2e-10 ***
## lytes1 0.9133 0.0384 23.80 < 2e-16 ***
## procedure 0.3482 0.0283 12.32 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(Age) 6.07 7.02 6621.6 < 2e-16 ***
## s(T) 4.14 5.11 468.9 < 2e-16 ***
## s(Nmonth) 3.69 10.00 57.6 2.4e-14 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.0631 Deviance explained = 8.13%
## -REML = 63575 Scale est. = 1 n = 196099
Logit
## [1] "AIC 127713.4842"
## [1] "RsqrMcFadden 0.0767"
## [1] "Residual deviance 127665.4842"
## [1] "Area under the ROC curve"
## [1] 0.70728
GAM Logit
## [1] "AIC 127089.5396"
## [1] "RsqrMcFadden 0.0813"
## [1] "Residual deviance 127035.3216"
## [1] "Area under the ROC curve"
## [1] 0.71062
Contraste sobre la reducción en la deviance
## Analysis of Deviance Table
##
## Model 1: Yvar ~ Age + Pub.Hol + chf + pulmcirc + chrnlung + renlfail +
## lymph + mets + coag + lytes + procedure + Trend + as.factor(Season)
## Model 2: Yvar ~ s(Age) + Pub.Hol + chf + pulmcirc + chrnlung + renlfail +
## lymph + mets + coag + lytes + procedure + s(T) + s(Nmonth,
## k = 12, bs = "cc")
## Resid. Df Resid. Dev Df Deviance Pr(>Chi)
## 1 196075 127665
## 2 196074 127035 0.904 630 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
El paquete mgcv proporciona un contraste sobre la calidad de los suavizados
## null device
## 1
##
## Method: REML Optimizer: outer newton
## full convergence after 5 iterations.
## Gradient range [-3.2083e-08,0.0010527]
## (score 63575 & scale 1).
## Hessian positive definite, eigenvalue range [0.43901,2.5099].
## Model rank = 39 / 39
##
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
##
## k' edf k-index p-value
## s(Age) 9.00 6.07 0.93 0.28
## s(T) 9.00 4.14 0.96 0.94
## s(Nmonth) 10.00 3.69 0.94 0.53
## null device
## 1
## (Intercept) Pub.Hol1 chf1 pulmcirc1 chrnlung1 renlfail1
## 0.088442 1.046001 2.317976 1.703086 1.179204 1.378106
## lymph1 mets1 coag1 lytes1 procedure
## 1.392714 5.024117 1.514758 2.492439 1.416518
## (Intercept) Age Pub.Hol1 chf1
## 0.0025 1.0552 1.0446 2.3624
## pulmcirc1 chrnlung1 renlfail1 lymph1
## 1.6828 1.1408 1.3694 1.3330
## mets1 coag1 lytes1 procedure
## 4.8804 1.4943 2.5325 1.3543
## Trend as.factor(Season)2 as.factor(Season)3 as.factor(Season)4
## 0.9967 0.9236 0.9430 0.8968
## as.factor(Season)5 as.factor(Season)6 as.factor(Season)7 as.factor(Season)8
## 0.8349 0.8483 0.8032 0.8625
## as.factor(Season)9 as.factor(Season)10 as.factor(Season)11 as.factor(Season)12
## 0.8199 0.8853 0.8808 0.9483
##
## Family: binomial
## Link function: logit
##
## Formula:
## Yvar ~ s(Age) + Pub.Hol + chf + pulmcirc + chrnlung + renlfail +
## lymph + mets + coag + lytes + procedure + s(T, by = proceditromF) +
## s(Nmonth, k = 8, bs = "cc")
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.4251 0.0108 -225.18 < 2e-16 ***
## Pub.Hol1 0.0451 0.0164 2.74 0.0061 **
## chf1 0.8405 0.0259 32.49 < 2e-16 ***
## pulmcirc1 0.5320 0.0465 11.45 < 2e-16 ***
## chrnlung1 0.1649 0.0241 6.84 8.0e-12 ***
## renlfail1 0.3203 0.0280 11.43 < 2e-16 ***
## lymph1 0.3312 0.1148 2.88 0.0039 **
## mets1 1.6146 0.0461 35.00 < 2e-16 ***
## coag1 0.4153 0.0676 6.14 8.1e-10 ***
## lytes1 0.9129 0.0384 23.79 < 2e-16 ***
## procedure 0.3654 0.0422 8.66 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(Age) 6.07 7.02 6622.0 < 2e-16 ***
## s(T):proceditromF0 3.55 4.39 435.7 < 2e-16 ***
## s(T):proceditromF1 3.82 4.73 41.3 8.0e-08 ***
## s(Nmonth) 3.31 6.00 56.8 3.2e-14 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.0631 Deviance explained = 8.14%
## -REML = 63576 Scale est. = 1 n = 196099
En nuestro modelo aditivo generalizado incorporamos los efectos aleatorios a nivel hospital (en la constante y en la variable “procedimiento”) estimando una nueva función de suavizado (de acuerdo con una distribución Normal) para la interacción entre la variable que identifica al hospital y la dummy de procedimiento.
##
## Family: binomial
## Link function: logit
##
## Formula:
## Yvar ~ s(Age) + Pub.Hol + chf + pulmcirc + chrnlung + renlfail +
## lymph + mets + coag + lytes + procedure + s(T, by = proceditromF) +
## s(Nmonth, k = 8, bs = "cc") + s(chos2, bs = "re", by = proceditromF)
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.4229 0.0396 -61.13 < 2e-16 ***
## Pub.Hol1 0.0456 0.0165 2.76 0.0058 **
## chf1 0.8332 0.0261 31.89 < 2e-16 ***
## pulmcirc1 0.5255 0.0468 11.22 < 2e-16 ***
## chrnlung1 0.1622 0.0242 6.69 2.2e-11 ***
## renlfail1 0.3211 0.0282 11.40 < 2e-16 ***
## lymph1 0.3262 0.1152 2.83 0.0046 **
## mets1 1.6293 0.0464 35.09 < 2e-16 ***
## coag1 0.4288 0.0679 6.31 2.8e-10 ***
## lytes1 0.8866 0.0388 22.88 < 2e-16 ***
## procedure 0.3881 0.0820 4.73 2.2e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(Age) 6.16 7.13 6611.8 < 2e-16 ***
## s(T):proceditromF0 3.05 3.79 433.9 < 2e-16 ***
## s(T):proceditromF1 3.98 4.91 41.5 8.0e-08 ***
## s(Nmonth) 3.22 6.00 54.0 2.2e-13 ***
## s(chos2):proceditromF0 34.43 36.00 890.4 < 2e-16 ***
## s(chos2):proceditromF1 25.23 36.00 134.2 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.0704 Deviance explained = 8.92%
## -REML = 2.7661e+05 Scale est. = 1 n = 196099
GAM Logit
## [1] "AIC 127089.5396"
## [1] "RsqrMcFadden 0.0813"
## [1] "Residual deviance 127035.3216"
## [1] "Area under the ROC curve"
## [1] 0.70728
Mixed GAM Logit
## [1] "AIC 126129.441"
## [1] "RsqrMcFadden 0.0892"
## [1] "Residual deviance 125943.3994"
## [1] "Area under the ROC curve"
## [1] 0.7188
Contraste sobre la reducción en la deviance
## Analysis of Deviance Table
##
## Model 1: Yvar ~ s(Age) + Pub.Hol + chf + pulmcirc + chrnlung + renlfail +
## lymph + mets + coag + lytes + procedure + s(T) + s(Nmonth,
## k = 12, bs = "cc")
## Model 2: Yvar ~ s(Age) + Pub.Hol + chf + pulmcirc + chrnlung + renlfail +
## lymph + mets + coag + lytes + procedure + s(T, by = proceditromF) +
## s(Nmonth, k = 8, bs = "cc") + s(chos2, bs = "re", by = proceditromF)
## Resid. Df Resid. Dev Df Deviance Pr(>Chi)
## 1 196069 127035
## 2 195995 125943 73.6 1092 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
El paquete mgcv dispone de una función “bam” para el ajuste de modelos GAM con grandes volúmenes de datos. La ventaja de bam es que optimiza el uso de la memoria en la computación utilizando menos tiempo para la convergencia de los modelos. bam también puede trabajar con clústers previamente configurados en el paquete parallet.